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Abstract — This  paper  proposes  a  geometry-based  statistical 
model  for  multiple-input  multiple-output  shallow  water  acous¬ 
tic  multipath  fading  channels.  From  the  reference  model,  the 
corresponding  space-time-frequency  correlation  function  is  de¬ 
rived.  Finally,  the  derived  spatial  and  temporal  correlations  are 
compared  with  the  empirically  obtained  channel  statistics  and 
close  agreement  is  observed1. 

I.  Introduction 

Shallow  water  acoustic  (SWA)  communications  find  ap¬ 
plications  in  oceanographic  data  collection,  tactical  surveil¬ 
lance,  and  offshore  exploration.  To  make  these  applications 
feasible,  there  is  a  need  for  wireless  underwater  acoustic 
communications  among  deployed  sensors  and  autonomous 
underwater  vehicles.  Statistical  characterization  of  acoustic 
communication  channels  is  necessary  in  order  to  assess  system 
performance  and  improve  the  quality  of  a  system  [1], 

The  SWA  channel  is  one  of  the  most  challenging  commu¬ 
nication  channels  because  it  suffers  path-dependent  Doppler 
and  angle  spreading,  which  results  in  a  time-varying  wideband 
channel  impulse  response  with  long  delay  spreads  [2].  While 
several  deterministic  and  statistical  simulation  models  of  SWA 
channel  have  been  proposed  [3]-[6],  a  statistical  framework 
for  SWA  channels  is  still  an  open  problem.  In  particular,  the 
space-time-frequency  correlation  function  is  a  necessary  tool 
for  a  proper  design  and  analysis  of  multiple-input  multiple- 
output  (MIMO)  SWA  communication  systems.  By  taking  into 
account  only  macro-scattering  effects,  Abdi  and  Guo  [7]  were 
the  first  to  derive  the  space-frequency  correlation  function  for 
time-invariant  isovelocity  SWA  channels. 

In  contrast,  this  paper  models  a  time-varying  SWA  channel 
by  taking  into  account  both  macro-  and  micro-scattering 
effects.  We  first  introduce  a  new  geometry-based  reference 
model  for  wideband  MIMO  SWA  multipath  fading  chan¬ 
nels.  The  proposed  model  characterizes  the  sound  propaga¬ 
tion  in  shallow  water  isovelocity  environments  by  combining 
the  deterministic  ray-tracing  theory  with  statistical  methods 
needed  to  characterize  random  components  of  the  propagation 
medium.  From  the  reference  model,  the  closed-form  space- 
time-frequency  correlation  function  (stf-cf)  for  a  shallow  water 
isovelocity  environment  is  derived.  The  closed-form  stf-cf 
can  be  useful  for  estimation  of  physical  parameters  of  the 
channel,  such  as  angle  spread,  mean  angles  of  departure  and 
arrival,  etc.  This  allows  system  engineers  to  understand  how 
these  channel  parameters  affect  the  correlation,  which  in  turn 
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provides  useful  guidelines  for  the  system  design.  To  illustrate 
the  utility  of  the  proposed  model,  we  compare  the  temporal 
and  spatial  correlation  functions  with  those  obtained  from  the 
measurements  reported  in  [8],  [9], 

The  remainder  of  the  paper  is  organized  as  follows.  Sec¬ 
tion  II  introduces  the  geometry-based  statistical  model  for 
wideband  MIMO  SWA  multipath  fading  channels.  Section  III 
derives  the  stf-cf  for  a  shallow  water  isovelocity  environment. 
Section  IV  compares  the  analytical  and  empirical  results 
for  the  temporal  and  spatial  correlations.  Finally,  Section  V 
provides  some  concluding  remarks. 

II.  Geometry-based  Statistical  Model  for 
Wideband  MIMO  SWA  Channels 

This  paper  considers  a  MIMO  communication  system  with 
Lt  transmit  and  Lr  receive  transducers.  The  propagation 
occurs  in  shallow  water  environments  with  a  constant  sound 
speed  and  is  characterized  by  two-dimensional  (2-D)  wide 
sense  stationary  uncorrelated  scattering  (WSSUS)  with  either 
line-of-sight  or  non-line-of-sight  conditions  between  the  trans¬ 
mitter  (Tx)  and  receiver  (Rx).  The  MIMO  SWA  channel  can 
be  described  by  an  Lr  x  Lt  matrix  H(f,  r)  =  \h-ij  (t:  t)] Lr x  Lt 
of  the  input-delay  spread  functions. 

Fig.  1  shows  a  SWA  channel  with  Lt  =  Lr  =  2  transducer 
elements.  This  elementary  2x2  transducer  configuration 
will  be  used  later  to  construct  uniform  linear  arrays  with  an 
arbitrary  number  of  transducer  elements.  The  SWA  channel 
is  modeled  as  a  2-D  waveguide  bounded  from  the  top  and 
bottom.  The  surface  and  bottom  boundaries  reflect  an  acoustic 
signal,  which  results  in  multiple  eigenrays  travelling  between 
the  Tx  and  Rx,  as  shown  in  Fig.  1.  The  roughness  of  sea 
surface  and  sea  bottom  is  characterized  by  S  and  B  macro- 
scatterers,  respectively.  At  any  time  instance  t,  the  Rx  receives 
2 S  downward  arriving  eigenrays,  each  one  having  different 
number  of  s  surface  and  b  bottom  reflections,  where  1  <  s  < 
S,  and  s  —  1  <  b  <  s.  Similarly,  there  are  2 B  upward  arriving 
eigenrays  with  b  bottom  and  s  surface  reflections,  where  1  < 
b  <  B  and  6—1  <  s  <  b.  Note  that  exact  positions  of  scatterers 
depend  on  the  surface  and  bottom  characteristics  and  may  vary 
from  one  location  to  another  and  from  one  time  instance  to 
another.  On  the  other  hand,  the  average  of  scattering  positions 
over  different  locations  (or  time  instances)  depends  only  on  the 
waveguide  geometry  and  the  number  of  eigenrays,  and  can  be 
computed.  To  implement  this  idea,  each  eigenray  is  modeled 
as  an  average  of  rays  scattered  from  Ns &  possible  positions 
of  surface  scatterers  and  Mhs  possible  positions  of  bottom 
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Fig.  1.  The  geometry-based  model  for  a  MIMO  wideband  SWA  channel 
with  Lt  =  Lr  =  2  transducer  elements. 


scatterers  (i.e.,  micro-scatterers),  where  the  possible  positions 
of  scatterers  are  clustered  around  the  averaged  positions  of 
scatterers,  as  shown  in  Fig.  1.  The  horizontal  spacing  between 
the  Tx  and  Rx  is  denoted  by  R.  The  water  depth  is  denoted  by 
h,  while  the  depths  of  Tx  and  Rx  are  denoted  by  hr  and  h,R, 
respectively.  As  we  are  interested  in  medium  and  long  range 
shallow  water  communications,  we  assume  that  the  depths  h , 
hr,  and  hp>  are  much  smaller  than  the  distance  R.  The  spacing 
between  two  adjacent  transducer  elements  at  the  Tx  and  Rx 
is  dr  and  da,  respectively.  It  is  assumed  that  dr  and  <Ir  are 
much  smaller  than  the  depths  h,  hr,  and  Hr.  The  symbols 
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eTsbn  an^  eTbsm 


denote  distances 


and  T. 


<J>)  q(B) 


Jbm 


respectively,  where  and  sj^  denote  the  nth  and  mth 
micro-scatterers  at  the  surface  and  bottom,  respectively,  for 
1  <  n  <  Ntb  and  1  <  to  <  M^t,  as  shown  in  Fig.  1. 
Similarly,  the  symbols  e. 
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Fig.  2.  The  detailed  geometry  of  the  DPP,  single-bounced  surface,  and  single- 
bounced  bottom  eigenrays  scattered  from  the  nth  and  mth  micro-scatterer, 
respectively. 


single-bounced  surface  and  single-bounced  bottom  eigenrays 
(i.e.,  s  =  1,  b  =  0,  s  =  0,  b  =  1)  scattered  from  the  nth  and 
mth  micro-scatterers,  respectively.  Angles  6t  and  Or  in  Fig.  2 
describe  the  orientations  of  Tx  and  Rx  transducer  arrays  in  the 
x-z  plane,  respectively,  relative  to  the  tr-axis.  The  Tx  and  Rx 
are  moving  with  speeds  vt  and  vr  in  directions  described 
by  angles  77-  and  7 r  in  the  x-z  plane  (relative  to  the  x- 
axis),  respectively.  The  symbols  cvf/,bn  anc^  aTbsm  are  t^le 
angles  of  departure  (AoD)  of  eigenrays  that  start  from  rr!,PI 

(  ( R') 

and  impinge  on  the  scatterers  Ssn  and  Sbm  ,  respectively, 
whereas  Oi/^sbn  and  o/rI  are  the  angles  of  arrival  (AoA)  of 
the  eigenrays  scattered  from  SXl  and  and  arriving  at 

R^x  \  respectively.  Finally,  the  symbols  aRdpp  and  D  denote 
the  AoA  of  DPP  ray  and  the  distance  Or-Or ,  respectively. 

Observe  from  the  geometrical  model  in  Fig.  1  that  prop¬ 
agation  through  SWA  channel  can  be  characterized  as  a 
superposition  of  surface-bottom  bounced  eigenrays  and  the 
DPP  eigenray  between  the  Tx  and  Rx.  The  multiple  bounced 
rays  can  be  grouped  into  the  upward  arriving  (UA)  rays  (i.e., 
the  last  reflection  is  from  the  bottom)  and  the  downward 
arriving  (DA)  rays  (i.e.,  the  last  reflection  is  from  the  surface). 
Then,  the  input  delay-spread  function  of  the  pressure  link  Tx  - 
Rx  can  be  written  as  a  superposition  of  the  DPP,  UA,  and 
DA  rays,  viz. 

hpq(t,T)  =  hpqpp(t,  r)  +  hp^(t,  r)  +  hpqA(t,  r).  (1) 

The  DPP  component  of  the  input  delay-spread  function  is 


hpq  (*> T)  V  K+1  VGp  (&Rddp  7r)  Gq  {oijidpp^ 

x  LS{D)LA{D)5(T-T%£\tj),  (2) 

where  Gp(-)  and  Gq ( ■ )  denote  the  radiation  patterns  of  the  pth 
transmit  and  qth  receive  transducer  element,  respectively,  K  is 
the  Rice  factor  (ratio  of  DPP  to  scatter  received  power),  Lr(D) 
and  LA(D)  denote  the  propagation  loss  due  to  spherical 
spreading  and  absorption,  respectively,  and  ()(■)  denotes  the 
Dirac  impulse  function.  The  DPP  time  delay  x//'^  (t)  denotes 
the  travel  time  of  the  direct  ray  between  Tx ^  and  R^\  i.e., 

Tdlf  (t)  =  +  l^tcos(aRdpp  +  tt  -  7 T) 

+  —  tcos(aRdpp  -  7ij),  (3) 

where  c  is  the  speed  of  sound.  We  assume  that  the  Tx  has 
omnidirectional  transducer  elements  and  therefore  produces  a 
spherical  wavefront  in  an  isovelocity  medium.  The  propagation 
loss  caused  by  the  spherical  spreading  can  be  written  as  [10] 
Lr(D)  =  1/D.  Furthermore,  when  sound  propagates  in  the 
ocean,  part  of  the  acoustic  energy  is  continuously  transformed 
into  heat.  This  absorption  is  primarily  a  result  of  relaxation 
processes  in  seawater.  The  absorption  loss  can  be  written  as 
[10] 


of  reference.  Fig.  2  details  the  geometry  of  the  deterministic 
(direct)  propagation  paths  (DPP)  as  well  as  the  geometry  of 


LA(D)  =  lQ-^/20000^ 


(4) 


where 


P  =  8.68  •  103  (5) 

x  +  BJ^{1  -  6'54 ' 10_4  '  P))  [dB/km], 

where  A  =  2.34  •  lO"6,  B  =  3.38  •  lO"6,  S  is  salinity  (in 
%o),  P  is  hydrostatic  pressure  [kg/cm2],  fc  is  the  carrier  fre¬ 
quency  [kHz]  and  fr  =  21.9  •  io6_152°/(T+273)  is  relaxation 
frequency  [kHz],  with  T  as  the  temperature  [°C]. 

The  upward  and  downward  arriving  components  of  the  input 
delay-spread  function  are,  respectively. 


hpaA{t,r)  = 
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n= 1 
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where  £bsm  >  0  and  £sbn  >  0  denote  the  amplitudes,  <j)sbn  and 
4>bsm  denote  the  phases  of  multipath  components,  Ls(Dsb) 
and  Ls(Dbs)  denote  the  propagation  loss  due  to  spherical 
spreading,  and  LA(Dsb)  and  LA(Dbs)  denote  the  propagation 
loss  due  to  absorption,  where  Dsb  and  Dbs  denote  the  total 
distances  travelled  by  UA  and  DA  rays  with  s  surface  and 
b  bottom  reflections.  These  distances  are  obtained  using  the 
method  of  images  and  can  be  written  as  [10] 


Dsb  =  ^R2  +  (2 bh  +  hR-  (-1  J°~b)hT)2,  (8) 
Dbs  =  ^R2  +  (2 sh  -hR  +  (-1  )(b-s)hT)2.  (9) 

The  parameters  rjs  and  tjb  in  (6)  and  (7),  respectively,  specify 
how  much  the  UA  and  DA  rays  contribute  in  the  total  power, 
i.e.,  these  parameters  satisfy  ijs  +  tjb  =  1.  The  parameters 
LsiPsh)  and  Lsjdbs)  denote  the  impedance  mismatch  be¬ 
tween  the  seawater  and  seabed  and  can  be  written  as  [10] 


Lb(0) 


m  cos  9  —  \J  n2  —  sin2  9 
m  cos  9  +  \J n?  —  sin2  9 


(10) 


where  m  =  pi/p,  n  =  c/cl,  p  and  c  denote  the  density  and 
sound  speed  in  the  seawater,  respectively,  whereas  p\  and  c\ 
denote  the  density  and  sound  speed  in  the  seabed,  respectively. 
The  incidence  angles  9sb  and  9bs  are  obtained  as  follows: 


bs 


=  tan 
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2sh  —  +  (— l)!b  s) hr 


(12) 


Finally,  the  time  delays  ( t)  and  r|f (t)  are  the  travel 
times  of  UA  and  DA  rays  with  s  surface  and  b  bottom 
reflections,  respectively,  i.e.. 
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where  A Zsbn(t)  and  A Zbsm(t)  denote  the  vertical  displace¬ 
ments  of  surface  scatterers  due  to  the  surface  motion  and  the 


distances  ersb,  tRsb,  CTbs.  and  tRbs  are  defined  in  Fig.  1. 

It  is  assumed  that  the  AoDs  {o^sbn  and  a^5ln),  the 
AoAs  (a(fllhn  and  cr^sm),  and  the  vertical  displacements 
A Zsbn(t)  and  A Zbsrn(t)  are  random  variables.  Furthermore,  it 
is  assumed  that  the  phases  (f>sbn  and  (f>bsm  are  uniform  random 
variables  on  the  interval  [— 7r,  tt)  that  are  independent  from  the 
AoDs,  the  AoAs,  and  the  vertical  displacements  at  the  surface. 

(p)  (p) 

Note  that  the  AoDs  {otfsbn  and  cdfbsm)  are  dependent  on  the 
AoAs  (o'En  and  a^bsrn),  respectively.  Assuming  small  angle 
spreads,  dr  -C  mm(h,T,h  —  hr),  and  9r  <C  mm(h.R,h  — 
Hr),  we  can  approximate  the  DoAs  and  AoAs  in  Fig.  2 

( p )  _  ,  ( p )  _  ,  ( p )  _  ,  (p)  _  , 

as  “ Tshn  ~ '  aTs6n  ~  aT  sbm  aTbsm  ~  aTbsm  ~  aTbsm > 

aRsbn  ~  aRsbn  ~  aRsbn,  and  aRbsm  ~  aRbsm  ~  aRbsm, 

where  aTsbn ,  aTbsm ,  aRsbn ,  and  aRbsrn  are  shown  in  Fig.  2. 
Assuming  that  each  ray,  when  interacting  with  the  surface 
and  bottom,  has  equal  incident  and  reflecting  angles,  we  can 
observe  that  o?Tsbn  —  tt  —  cxRsbn  and  ctThsm  37t  cv-Rbsm- 
The  distances  e^bn,  e^sm,  e^sbn,  e^sm,  and  epq  can  be 
expressed  as  functions  of  the  random  angles  aRsbn  and  aRbsm 
and  the  angle  ocRdpp  as  follows 
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(19) 


where  parameters  p  and  q  take  values  from  the  sets  p  £ 
{1  and  q  £  {1  respectively.  The  deriva¬ 

tions  of  approximations  in  (15)-(19)  are  omitted  for  brevity. 
To  simplify  further  analysis,  we  normalize  all  time  delays 
with  respect  to  the  direct  arrival  time.  This  normalization  is 
reflected  in  the  distances  epq,  Dsb,  and  Dbs.  Furthermore, 


since  the  depths  h,  hr,  and  h R  are  much  smaller  than  the 
distance  R,  the  distances  Dsb  and  Dbs  can  be  approximated 
using  y'l  +  x  «  1  +  x/2,  for  small  x.  Hence,  the  distances 
epq,  Dsb,  and  Dbs  can  be  approximated  as 


epq  «  0.5(Lt  +  1  —  2 p)dr  cos(a_R dpp  —  Qt) 

-0.5 (Lr  +  1  -  2 q)dR  cos{aRdpp  -  dR),  (20) 

Dsb  «  hb2h2  +  bhhR-{-l)^bhhT  +  {s  -  b)hThR\  (21) 
K 

Dbs  ~  -kl^h? -  shhR  +  {-l)^s) shhT  +  (b  -  s)hThR\.{ 22) 

H 


Finally,  we  note  that  ersb,  £Rsb,  tTbs*  and  eRbs  can  be  written 
as  functions  of  the  mean  AoA  angles  pRsb  and  pRbs,  he., 
tTsb  =  hT/ sin  nRsb,  eRsb  =  hR/sm/j,Rsb,  eTbs  =  —(h — 
hr)/  sin  Hubs,  and  eRbs  =  -{h  -  hR)/sm/jRbs,  where  the 
mean  AoA  angles  pRsb  and  pRbs  are  depicted  in  Fig.  1. 

The  angles  of  arrival  aRsbn  and  aRbsm  are  modeled  using 
the  following  Gaussian  probability  density  functions  (pdfs) 
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for  7 r  <  aRbsm  <  2ir, 
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III.  Space-Time-Frequency  Correlation  Function 
of  the  Statistical  Model 


where  pRsb  and  pRbs  denote  the  the  mean  AoAs,  whereas 
aRsb  and  crRbs  denote  the  angle  spreads. 

Finally,  the  vertical  displacements  AZsbn[t)  and  AZbsm(t) 
are  modeled  as  zero-mean  Gaussian  random  processes  with 
stationary  and  independent  increments,  i.e.. 


f(AZsbn(t)) 

f(A-Zbsm(t)) 


exp 


(A  Zsbn(t))2 
(2^Ca  zsb) 


,  (25) 
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exp 


(A  Zbsm{t))2 
(2iCi  zbs) 


,(26) 


where  tC,\Zsb  and  t(j±Zb  denote  variances. 

To  simplify  further  analysis,  we  use  the  time-variant  transfer 
function  instead  of  the  input  delay-spread  function  and  we 
normalize  the  gain  patterns  of  the  transducer  elements  to  unity 
(i.e.,  we  assume  omnidirectional  array  elements),  although 
other  gain  patterns  can  be  accommodated  at  this  point.  The 
time-variant  transfer  function  is  the  Fourier  transform  of  the 
input  delay-spread  function  and  can  be  written  as 

Ppq(t,f)  =  FT{hpq{t,T)}  =  P°pp{t,f) 

+  PPA(t,f)  +  PPA(t,f)  (27) 

where  PpDqpp(t,f),  PPA(t,f),  and  PpDA(t,f)  are  DPP,  UA, 
and  DA  components  of  the  time-variant  transfer  function, 
respectively.  Using  (2)-(22),  the  DPP,  UA,  and  DA  components 
of  the  time- variant  transfer  function  in  (27)  can  be  written  as 


Ppqpp(t,f )  =  Fr  { hppp(t,T)}  =  ^ J^-1LS{D)LA{D) 

cos {aRdpp-'yT)-VRt  cos(aRdpp-'yR)] 


Assuming  a  2-D  WSSUS  isovelocity  shallow  water  en¬ 
vironment,  we  here  derive  the  stf-cf  of  the  geometry-based 
statistical  model.  The  normalized  stf-cf  between  two  time- 
variant  transfer  functions  Ppq(t,  f)  and  Ppq(t  +  At,  f  +  A/), 
is  defined  as 


Rpq,pq(At,Af) 


E[Ppq(t,f)*Ppq(t  +  At,f  +  Af)] 

vmutJWjmWJWi 


where  (  •  )*  denotes  complex  conjugate  operation,  E[  -  ]  is  the 
statistical  expectation  operator,  p,p  £  {1, . . . ,  Lt},  and  q,q  £ 
{1, . . . ,  Lr}.  Since  PpqA(t,  f)  and  PpqA{t ,  f)  are  independent 
zero-mean  random  processes,  (31)  can  be  simplified  to 


RpqM{At,Af)  =  RpPP(At,  A/)  (32) 

+  RupqApq{ At ,  A/)  +  RPApq(At,  A/) , 

where  RppF~(At,  A/),  RPA~~(At,  A/),  and  RPA~~(At,  A/) 

denote  the  normalized  stf-cfs  of  the  DPP,  UA,  and  DA  com¬ 
ponents,  respectively,  and  are  defined  as 
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333) 
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where  Qua  =  E[Ls(Dbs)LA(Dbs)LB{Obs)b],  Qdpp  = 
Ls{D)La{D),  and  QDA  =  E{Ls{Dsb)LA{Dsb)LB{esb)b]. 


By  substituting  (28)  into  (33),  the  expression  for  the  stf-cf 
of  DPP  component  becomes 

RpqPpPq(At’Af)  =  K 

KAi  cos(aRdpp—7T)—VRAt  cos(aRdpp-jR)] 
ei [(p-p)dr  cos(aRdpP-9T)-(q-q)d.R  cos(aRdpp-0R)] 


Substituting  (29)  and  (30)  into  (34)  and  (35),  respectively, 
noting  that  the  phases  (j>bsm  and  (j>sbn  are  independent  and 
uniformly  distributed  over  [— n,  7r),  and  assuming  that  Mbs 
1  and  Nsb  1,  the  stf-cfs  of  the  UA  and  DA  components 
can  be  written  as,  respectively, 
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For  small  angle  spreads,  the  AoAs  a#sb  and  are 

mainly  concentrated  around  the  mean  AoAs  and 
respectively.  Using  the  first-order  Taylor  expansion,  the  AoA 
angles  can  be  approximated  as  follows 

cos(a_RSb)  «  cos (uRsb)  -  sin(fj,Rsb)(aRsb  —  HRsb),  (39) 
sin(afisb)  «  sm(nRsb)  +  cos(nRsb)(aRsb  -  HRsb),  (40) 
sin(a_Rsb)2  «  sm(iiRSb)2 

+2sin(^fish)  cos(nRSb)(aRSb  —  HRsb ),  (41) 

sin(a_Rs6)_1  «  sin^Rsb)-1 

—  cos(fiRsb)(aRsb  -  HRsb)/ sm(nRsb)~2 ■  (42) 


The  similar  approximations  are  applied  to  the  AoAs  aRbs. 
Using  these  trigonometric  approximations  and  the  equality 
J  e^°‘x(2TTa2)~1^2e~x  R2a  ^dx  =  e_<T  “  /2  [7],  the  stf-cfs  of 
UA  and  DA  components  can  be  closely  approximated  as 
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IV.  Comparison  with  Measured  Data 
In  this  section,  we  compare  the  theoretical  results  in  Sec¬ 
tion  III  with  the  AUVFest07  measured  data  in  [8],  collected 
at  rough  sea.  Furthermore,  we  compare  the  derived  spatial 
correlation  function  with  the  spatial  correlation  functions  of 
the  exponential  model  in  [9],  the  particle -velocity  model  in  [7], 
and  the  measured  data  in  [9].  The  close  agreement  between 
the  analytical  and  measured  statistics  is  observed. 

The  channel  measurements  in  [8]  are  collected  at  fc  = 
17  kHz  and  the  speed  of  sound  was  c  =  1470  m/s.  The 
distance  between  the  Tx  and  Rx  was  R  =  5  km.  The  water, 
Tx,  and  Rx  depths  were  h  =  20  m,  hr  =  19  m,  and 
Hr  =  18  m,  respectively.  The  Tx  and  Rx  are  equipped  with 
vertically  oriented  transducer  elements  (i.e..  Or  =  Or  =  90°). 
It  is  assumed  that  the  Tx  and  Rx  are  relatively  stationary,  i.e., 
slightly  moving  with  waves  (i.e.,  7r  =  7 r  =  90°). 


Fig.  3.  The  magnitude  of  the  theoretical  and  measured  temporal  correlation 
functions  for  rough  sea. 


Fig.  3  compares  the  analytical  temporal  correlation  function 
in  (32)  with  the  measured  temporal  correlation  function  in  [8] 
for  the  rough  sea  environment.  The  analytical  curve  is  obtained 
with  the  estimated  parameters  K  =  0.38,  S  =  1,  B  =  1,  vr  = 
0.005  m/s,  vR  =  0.002  m/s,  Caz.=1,6=0  =  CAZa=1,b=i  = 
0.0014,  (Azb=1  a=0  =  CAZb=i  a=i  =  0.006,  fj,a~  1,6=0  =  168°, 
AU=i,6=i  =  171°,  Hb=i,s=o  =  205°,  Mfe=M=i  =  210°, 
crs=l,fc=0  =  2.55°,  £78=1,6=1  =  3.1°,  £76=1)S=0  =  10°, 
<Jb- i)S=i  =  12°,  r]s  =  0.2,  and  tjr  =  1— 77 s-  The  parameters  S 
and  B  are  visually  estimated  from  Fig.  2d  in  [8].  The  rest  of 
the  parameters  are  estimated  using  the  maximum  likelihood 


method  described  in  [12],  The  parameters  R,  h,  Or,  Hr, 
7 t,  'Jr,  Or,  Or,  and  ctRdpp  «  tt,  are  selected  to  match  the 
measurement  conditions  as  described  above.  Fig.  4  compares 


Fig.  4.  The  magnitude  of  the  theoretical  and  measured  spatial  correlation 
functions  for  rough  sea. 

the  analytical  spatial  correlation  function  in  (32)  with  the 
measured  spatial  correlation  function  in  [8]  for  the  rough  sea 
environment. The  analytical  curve  in  Fig.  4  is  obtained  using 
the  same  parameters  as  in  Fig.  3.  Finally,  Fig.  5  compares  the 
spatial  correlation  in  (32)  with  the  analytical  spatial  correlation 
in  [7]  and  the  analytical  and  empirical  spatial  correlations  in 
[9].  The  ASCOTOl  channel  measurements  are  collected  at 


Antenna  Element  Number 

Fig.  5.  Comparison  of  the  analytical  spatial  correlations  in  (32),  [7],  and  [9] 
with  the  empirical  spatial  correlation  in  [9], 

the  carrier  frequency  1.2  kHz  and  the  speed  of  sound  was 
c  =  1440  m/s  (i.e.,  A=1.2  m)  [9].  The  distance  between  the  Tx 
and  Rx  was  li  =  10  km.  The  water,  Tx,  and  Rx  depths  were 
h  =  103  m,  Iit  =  50  m,  and  Iir  =  50  m,  respectively.  The  Rx 
was  equipped  with  33  vertically  oriented  transducer  elements 
(i.e.,  Or  =  90°)  with  L  =  Or  =  0.5  m  element  spacing, 
while  the  Rx  was  equipped  with  one  vertical  transducer  (i.e.. 
Or  =  90°).  The  vertical  correlation  is  calculated  with  respect 


to  eighth  element  from  the  bottom  of  the  3  3 -element  array.  The 
analytical  spatial  correlation  in  (32)  is  obtained  with  the  esti¬ 
mated  parameters  K  =  0.71,  S  =  B  =  1,  fis- 17=0  =  173.8°, 
/xs=i,b=i  =  173.4°,  /q,=i,s=o  =  183°,  pb= i>s=i  =  184.5°, 
<78=17=0  =  9.1°,  07=17=1  =  10.2°,  o7=iiS=o  =  5.59°, 
o'b=ijS=i  =  8.89°,  r/s  =  0.44,  and  r]B  =  1  —  775.  The 
parameters  R,  h,  hr,  Or,  Or,  Or,  and  a,Rdpp  «  7r,  are  selected 
to  match  the  measurement  conditions  described  above.  For  the 
reference.  Fig.  5  also  shows  the  spatial  correlation  functions 
of  the  exponential  model  in  [9],  i.e,  exp{— L2/(2A)2}  and  the 
particle-velocity  model  in  [7]  (with  parameters  Ab  =  0.56, 
/ib  =  183°,  =  173°,  <j, b  =  2.29°,  and  as  =  8°).  One  can 

observe  that  the  proposed  model  provides  a  closer  match  with 
the  experimental  correlation  than  the  models  in  [7]  and  [9]. 

The  close  agreement  between  the  analytical  and  measured 
statistics  in  Figs.  3-5  confirms  the  utility  of  the  proposed 
model. 

V.  Conclusions 

This  paper  proposed  the  geometry-based  statistical  model 
for  MIMO  SWA  fading  channels.  From  the  statistical  model, 
the  corresponding  space-time -frequency  correlation  function  is 
derived.  Finally,  the  derived  statistics  are  compared  with  the 
empirically  obtained  channel  statistics  and  the  close  agreement 
is  observed. 
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